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Abstract 

In  this  paper,  we  consider  the  issue  of  convergence  toward  entropy  solutions  for  high 
order  finite  volume  weighted  essentially  non-oscillatory  (WENO)  scheme  and  discontinuous 
Galerkin  (DG)  finite  element  method  approximating  scalar  nonconvex  conservation  laws.  Al¬ 
though  such  high  order  nonlinearly  stable  schemes  can  usually  converge  to  entropy  solutions 
of  convex  conservation  laws,  convergence  may  fail  for  certain  nonconvex  conservation  laws. 
We  perform  a  detailed  study  to  demonstrate  such  convergence  issues  for  a  few  representative 
examples,  and  suggest  a  modification  of  the  high  order  schemes  based  either  on  first  order 
monotone  schemes  or  a  second  order  entropic  projection  [1]  to  achieve  convergence  toward 
entropy  solutions  while  maintaining  high  order  accuracy  in  smooth  regions. 
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1  Introduction 


In  this  paper,  we  consider  the  Cauchy  problem  of  one  dimensional  scalar  conservation  laws: 

f  ut  +  f(u)x  =  0,  in  R  x  [0, T\, 

\  u(x,  0)  =  u0(x),  in  R.  '  '  ’ 

The  unique  entropy  solution  of  (1.1)  satisfies 


U(u)t  +  F(u)x  <  0 


(1.2) 


in  the  distribution  sense,  for  any  convex  entropy  function  U (■ u )  and  its  corresponding  entropy 
flux  function  F(u)  satisfying  F'(u)  =  U'(u)f'(u). 

If  (1.1)  is  a  convex  (or  concave)  conservation  law,  that  is,  if  f"(u)  does  not  change  sign, 
we  would  only  need  the  entropy  condition  (1.2)  to  be  satisfied  for  one  strictly  convex  entropy 
in  order  to  determine  the  unique  entropy  solution  of  (1.1).  When  (1.1)  is  non-convex,  namely 
when  f"{u )  changes  sign,  the  solution  structure  of  (1.1)  is  much  more  complicated.  This  can 
already  be  seen  from  the  simple  Riemann  problem  with  the  initial  condition 


u0(x) 


Ui,  for  x  <  0, 

Ur,  for  X  >  0. 


(1.3) 


The  discontinuity  at  the  origin  generates  a  single  shock  or  a  rarefaction  wave  for  a  convex 
conservation  law,  while  for  nonconvex  conservation  laws,  much  more  complicated  solution 
structures,  sometimes  called  a  compound  wave,  might  be  generated.  The  compound  wave 
involves  a  sequence  of  shocks  and  rarefaction  waves  and  is  much  harder  to  be  resolved 
numerically.  The  entropy  solution  of  the  nonconvex  Riemann  problem  can  be  determined 
from  a  convex-hull  construction  [8]. 

Numerical  methods  for  nonlinear  hyperbolic  equations  are  expected  to  be  able  to  ef¬ 
ficiently  capture  the  unique  entropy  solution.  It  is  well  known  that  first  order  monotone 
schemes  converge  to  entropy  solutions  of  both  convex  and  nonconvex  conservation  laws  [3], 
but  with  a  relatively  slow  convergence  rate.  For  high  order  numerical  methods,  for  example 
the  high  order  finite  difference  or  finite  volume  WENO  schemes  [9,  6,  13]  and  discontinuous 
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Galerkin  method  [2],  however,  very  few  theoretical  results  on  convergence  are  available  for 
discontinuous  solutions. 

The  high  order  finite  volume  WENO  scheme,  about  which  we  provide  a  brief  review  in 
section  2.1,  is  one  of  the  successful  high  order  numerical  methods  in  approximating  (1.1). 
Based  on  the  idea  of  adaptive  stencils  in  the  reconstruction  procedure,  WENO  schemes  would 
automatically  achieve  high  order  accuracy  in  smooth  regions  of  the  solution  and  an  essentially 
non-oscillatory  resolution  at  the  discontinuities.  Despite  its  good  numerical  performance  in 
many  applications,  the  theoretical  proof  for  its  convergence  towards  the  entropy  solution 
is  very  difficult.  Typically,  one  would  need  a  cell  entropy  inequality  [10],  or  a  wavewise 
entropy  inequality  [17],  together  with  the  total  variation  bounded  (TVB)  property  of  the 
numerical  solutions  to  imply  such  convergence.  There  are  counter  examples  for  nonconvex 
conservation  laws  [18] ,  such  that  high  order  Godunov  schemes  would  not  be  able  to  converge 
to  the  entropy  solution  under  the  usual  time  step  governed  by  the  CFL  condition.  For  the 
Godunov  schemes  under  large  time  steps,  convergence  is  available  with  up  to  third  order 
reconstruction  operators  [11].  Though  this  convergence  result  holds  for  both  convex  and 
nonconvex  conservation  laws  in  multiple  space  dimensions,  the  implementation  of  the  scheme 
for  nonconvex  conservation  laws  and  for  multiple  space  dimensions  is  highly  non-trivial. 

The  discontinuous  Galerkin  (DG)  method  is  another  class  of  high  order  numerical  meth¬ 
ods  in  approximating  (1.1),  about  which  we  provide  a  brief  review  in  section  2.2.  ft  is  a  finite 
element  method,  but  uses  discontinuous,  typically  piecewise  polynomial  functions  to  form 
the  solution  and  test  spaces,  and  adopts  finite  volume  techniques  including  numerical  fluxes 
and  limiters.  It  can  be  proved  that  the  semi-discrete  DG  method  satisfies  a  cell  entropy  in¬ 
equality  for  the  square  entropy  [5],  which  implies  convergence  to  the  unique  entropy  solution 
for  one  dimensional  scalar  convex  conservation  laws  if  it  converges.  For  nonconvex  conserva¬ 
tion  laws,  an  entropy  inequality  for  a  single  entropy  is  not  enough  to  imply  convergence  to 
the  entropy  solution.  Examples  for  which  numerical  solutions  of  the  DG  method  fail  to  con¬ 
verge  to  the  correct  entropy  solution  for  nonconvex  conservation  laws  will  be  demonstrated 
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in  section  3. 


In  this  paper,  we  propose  a  first  order  monotone  modification  in  the  framework  of  high 
order  finite  volume  WENO  scheme  and  the  discontinuous  Galerkin  method  to  enforce  con¬ 
vergence  towards  the  entropy  solution.  A  discontinuity  indicator  is  introduced  to  avoid  the 
degeneracy  of  high  order  accuracy  in  smooth  regions.  Based  on  the  idea  of  a  MUSCL  type 
scheme  satisfying  all  the  entropy  inequalities  in  [1],  a  second  order  modification  with  an 
entropic  projection  is  also  proposed  for  high  order  numerical  schemes. 

This  paper  is  organized  as  follows.  Section  2  gives  a  brief  review  on  high  order  finite  vol¬ 
ume  WENO  and  discontinuous  Galerkin  methods.  Section  3  demonstrates  the  performance 
of  different  numerical  methods  for  nonconvex  conservation  laws  through  several  examples. 
A  first  order  monotone  modification,  as  well  as  a  second  order  modification  with  an  entropic 
projection  are  proposed  for  high  order  numerical  methods  in  section  4  and  section  5  re¬ 
spectively.  Numerical  examples  are  shown  to  demonstrate  the  quality  of  proposed  schemes. 
Concluding  remarks  are  given  in  section  6. 

2  High  order  numerical  methods 

In  this  section,  we  briefly  review  two  classes  of  high  order  numerical  schemes,  i.e.  the  high 
order  finite  volume  WENO  scheme  and  the  discontinuous  Galerkin  method  in  approximating 
(1.1).  This  forms  the  basis  for  our  discussion  in  later  sections  on  high  order  numerical 
methods.  Throughout  the  paper,  we  will  always  use  the  following  notation  for  the  meshes 
unless  otherwise  specified.  The  computational  interval  is  divided  into  N  subintervals,  each  of 
which  is  denoted  as  a  cell  Ij  =  [Xj_i,Xj+i\  with  the  cell  center  Xj  =  +Xj+i)  and  the 

cell  length  A  Xj  =  (xJ+i  —  Xj_i).  We  also  denote  the  maximum  cell  size  as  Ax  =  max,,-  A  Xj. 
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2.1  The  finite  volume  WENO  scheme 


One  of  the  successful  high  order  numerical  methods  for  (1.1)  is  the  class  of  finite  volume 
WENO  method  [9,  6,  13].  It  approximates  the  integral  version  of  (1.1)  by 


du3  _ _ L  (f  1 

dt  AXj^j+* 


(2.1) 


where  Uj  is  the  approximation  to  the  cell  average  of  the  solution  over  the  cell  L 


ui  —  — 

J  Ax 


u{x ,  t)dx, 


j  Jii 


and 


fj+i  =  /K+n<a) 


(2.2) 


j+r  j+3 

is  a  monotone  numerical  flux  (non-decreasing  in  the  first  argument  and  non- increasing  in 
the  second  argument)  at  the  cell  boundary  x]+\ .  Examples  of  monotone  fluxes  include  the 
Godunov  flux, 

?G(„  u\  _  )  ^  3) 


mi n  /(0,  if  a  <  b, 

fG(a,b)={  *e[aMfU,  ^ 

'  j  max  /(£),  if  b  <  a, 
«e[M 


which  is  the  least  dissipative  among  all  monotone  fluxes,  and  the  Lax- Friedrichs  flux, 


f  (a,  b)  =  ~{f{a)  +f(b)  -  a(b-a)),  a  =  max|/'(u)|. 


(2.4) 


In  (2.2),  u~+i  and  u++i  are  the  reconstructed  values  of  the  solution  from  neighboring  cell  av¬ 
erages  to  the  left  and  right  of  the  cell  interface  x]+i  ■  Examples  of  high  order  reconstructions 
include  second  order  MUSCL  reconstruction  and  fifth  order  WENO  reconstruction. 


•  Second  order  MUSCL  reconstruction  [15]. 

The  reconstructed  point  values  for  the  second  order  MLISCL  scheme  are  obtained  by 

1  ,  _  1 

V*  =  ui  +  2S j)  V*  =  Uj+1  ~  2Sj+u  ( 2 -5) 

where  Sj  =  m{uj+ 1  —  Uj,Uj  —  Uj-i)  with  the  usual  dehnition  of  the  minmod  function 
m  [4] 

m(a,b )  =  -(sign(a)  +sign(6))  min(|a|,  \b\). 
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Fifth  order  WENO  reconstruction  [6,  13]. 

The  reconstructed  point  value  for  the  fifth  order  WENO  scheme  is  obtained  by 


“j+1/2  =  UlU)+l/2  +  ^Uj+l/2  +  ^3  u)+1/2,  (2-6) 


where  tv+j  ,,2  are  the  three  third  order  reconstructions  on  three  different  stencils  given 
by 

(D  _I-  _Z- 

'"7+1/2  —  3Mj_2  QUj~1  +  6  Uj ’ 

(2)  1  5  _  1 

u)U/2  =  -g%-l+  g%  +  3%+l, 

(3)  1_  5  1_ 

Mj+i/2  —  3%  +  g^H1  6Mj+2’ 

and  the  nonlinear  weights  u Jm  are  given  by 

&m  ~  _  Ti 

—  FT?  T" ;  ui  ~  ,  o\2  ’ 

2si=1vi  {£  +  Pi) 

with  the  linear  weights  7 1  given  by 


and  the  smoothness  indicators  given  by 

13  _  2  1  _  _  2 

Pi  =  ^2  (^'“2  “  2ui~l  +  +  4  (Mi-2  _  +  3wi) 

P2  =  ^  («i-i  -  2%  +  -Uj+r)2  +  ^  («J-1  -  %+i)2 

13  _  2  1  -  -  2 

/?s  =  ^2  ^  ~~  2“i+1  +  +  4  _  4“i+1  +  Uj+ 2’  ' 

Here  £  is  a  parameter  to  avoid  the  denominator  to  become  0  and  is  taken  as  10-6. 

The  reconstructed  point  value  u++ , ^  is  obtained  in  a  mirror  symmetric  fashion  with 

respect  to  37+1/2  as  that  for  uj+l,2.  For  more  details  on  the  reconstruction  procedure, 

we  refer  to  [6,  13]. 


The  high  order  Godunov  type  scheme  (see,  e.g.  [18,  11]  for  the  definition  of  this  scheme 
that  we  are  using)  is  a  fully  discretized  version  of  the  finite  volume  scheme  (2.1).  It  consists 
of  the  following  three  stages  to  evolve  the  cell  averages  from  un  at  n-tli  time  level  to  un+ 1 : 
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1.  Reconstruction:  obtain  a  high  order  piecewise  polynomial  reconstruction  un(x)  whose 
cell  averages  agree  with  the  given  cell  averages  un.  We  denote  this  reconstruction 
operator  as  Re(un).  Second  order  MUSCL  and  fifth  order  WENO  reconstructions  as 
described  above  are  typical  examples  of  this  reconstruction  operator. 

2.  Evolution:  evolve  un(x)  by  the  conservation  law  (1.1)  exactly  for  a  time  step  At,  to 
obtain  a  solution  un+1(x)  which  is  in  general  not  piecewise  polynomial  anymore.  We 
denote  this  evolution  operator  as  S^t(un). 

3.  Averaging:  average  the  function  un+l(x)  to  obtain  the  cell  averages  un+1  at  time  level 
n  +  1.  We  denote  this  averaging  operator  as  A(un+1). 

The  Godunov  type  scheme,  using  the  notations  introduced  above,  can  be  described  ab¬ 
stractly  as 

j  u°  =  A(u0)  ( 

\  un+1  =  AoSAtoRe(un),  n  =  0, 1,  •  •  •  1  ’ 

We  remark  that  the  exact  evolution  step  is  difficult  to  implement,  however  we  do  not 
need  all  the  information  of  this  exactly  evolved  solution  un+1(x),  but  only  its  cell  averages. 
Therefore,  it  is  usually  possible  to  implement  such  Godunov  type  schemes  in  an  efficient 
fashion.  We  refer  to  [11]  for  more  details. 

2.2  The  discontinuous  Galerkin  method 

Another  class  of  successful  high  order  numerical  methods  for  (1.1)  is  the  discontinuous 
Galerkin  method.  The  discontinuous  Galerkin  method  is  a  finite  element  method  using 
discontinuous  piecewise  polynomial  functions 

VL  =  {P  :  P\r,  e  Pk(IJ),  j  =  1,  2, ...,  N},  (2.8) 

as  the  solution  space  and  the  test  space,  where  Pk(Ij )  denotes  the  set  of  polynomials  of  degree 
<  k  on  cell  Iy  If  we  adopt  a  local  orthogonal  basis  over  Iy  say  {v^\x),  l  —  0,1, ...,  k },  then 
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the  numerical  solution  uAx(x,t )  G  V£x  can  be  written  as 


k 

uAx{x,t)  =  y (t)vij\x)  for  x  G  Ij 
1=0 


(2.9) 


with 

u^\t)  —  —  [  uAx{x,t)tp\x)dx,  l  =  0, 1,  2, k,  (2-10) 

ai  Jh 

where  a;  =  Jf  (y^\x))2dx  are  the  normalization  constants. 

A  semidiscrete  DG  method  for  solving  (1.1)  is  obtained  by  multiplying  (1.1)  by  a  test 
function  v(x)  G  V£x,  integrating  over  Ij,  and  integrating  by  parts.  Specifically,  the  evolution 
of  the  degrees  of  freedom  up  (t)  in  (2.9)  by  the  DG  method  is  achieved  by 


d  (,)  1 

~nuj  S 
dt  J  ai 


ix^))-^.vi\x)dx  +  fj+$vi(xj+$)  -  fi-^l\Xi- 


i)l  =0,  (2.11) 

0,1,-,  A:. 


Here,  fJ+i  is  a  monotone  numerical  flux  defined  by  (2.2)  with  and  v=+i  being  the  left 
and  right  limits  of  the  discontinuous  solution  uAx  at  the  cell  interface  1 . 

The  DG  method  is  a  linear  method,  hence  it  cannot  capture  discontinuities  without 
spurious  oscillations  when  the  order  of  accuracy  is  higher  than  one.  Typically,  nonlinear 
limiters  similar  to  those  used  in  high  resolution  finite  volume  schemes  are  used  to  suppress 
those  oscillations.  In  this  paper,  we  adopt  the  WENO  limiter  proposed  in  [12],  We  refer  to 
[12]  for  the  details  of  the  implementation  of  this  WENO  limiter. 

2.3  Total  variation  diminishing  (TVD)  Runge-Kutta  time  dis¬ 
cretization 

To  evolve  Uj(t)  in  (2.1)  for  the  finite  volume  WENO  scheme  or  ud\t)  in  (2.11)  for  the  DG 
method  in  time,  we  use  the  third  order  total  variation  diminishing  (TVD)  Runge-Kutta 
method  in  [14],  Explicitly,  for  the  system  of  ODEs  ut  =  L{u),  the  evolution  from  un  to  un+1 
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is  obtained  through 


«(1)  =  un  +  AtL(un,tn) 

u{2)  =  ^un  +  ^(u(1)  +  AtL(uw))  (2.12) 

un+1  =  \un  +  \{u^  +  AtL(u^)). 
o  o 

The  details  of  its  implementation  can  be  found  in  [14]. 


3  Examples  of  scalar  nonconvex  conservation  laws 

3.1  Examples  with  good  performance  of  high  order  schemes 

In  this  subsection,  we  will  give  two  examples  of  nonconvex  conservation  laws  for  which  the 
high  order  numerical  methods  perform  nicely  and  converge  to  the  entropy  solutions. 


The  Buckley-Leverett  Equation.  This  is  a  one-dimensional  model  problem  for  two- 
phase  fluid  flow  in  a  porous  medium  with  the  application  in  oil-reservoir  simulation. 
The  conservation  law  is  (1.1)  with  the  flux  function  given  by 


/(«)  = 


u 


v?  +  a(l  —  u )2 
and  the  initial  condition  given  by 


u0(x)  = 


1, 

0, 


for 

for 


where  a  =  -, 


x  <  0, 

x  >  0. 


(3.1) 


(3,2) 


Figure  3.1  shows  the  solution  u(x,t )  of  the  Buckley-Leverett  equation  (1.1)-(3.1)  with 
the  initial  condition  (3.2)  at  the  time  t  =  2.  The  exact  solution  is  obtained  by  the 
convex- hull  construction  [8].  We  observe  that  both  the  finite  volume  WENO  scheme 
and  the  discontinuous  Galerkin  method  with  a  WENO  limiter  approximate  the  exact 
solution  very  well  in  this  example.  We  have  computed  this  example  using  the  discon¬ 
tinuous  Galerkin  method  with  a  WENO  limiter  of  orders  of  accuracy  from  2  to  4  and 
have  obtained  similarly  good  solutions,  even  though  only  the  second  order  result  is 
plotted  in  Figure  3.1  to  save  space. 
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Figure  3.1:  Solid  lines:  the  exact  solution  of  (1. 1)- (3. 1)  with  the  initial  condition  (3.2)  at 
the  time  t  =  2;  solid  symbols:  the  numerical  solution  of  the  fifth  order  finite  volume  WENO 
scheme  (left)  and  the  second  order  discontinuous  Galerkin  method  with  WENO  limiters 
(right).  The  uniform  mesh  size  is  Ax  =  0.1. 


•  Another  example.  The  following  example  is  used  in  [7].  It  is  a  Riemann  problem  for 
the  scalar  one-dimensional  nonconvex  conservation  laws  (1.1)  with  the  flux  function 
given  by 


/(«)  =  < 

u{  1  —  u) 

4 

1  2  1 
-u - u  + 

l  2  2 

3 

16* 

.,  i 

11  u  <  -, 

K 

if  U>-, 

(3.3) 

and  the  initial  condition  given  by 

u0(x)  =  | 

for 

for 

x  <  0, 

x  >  0. 

(3.4) 

Figure  3.2  shows  the  solution  u(x,t)  for  the  Riemann  problem  of  (l.l)-(3.3)  with  the 
initial  condition  (3.4)  at  the  time  t  =  2.  The  exact  solution  is  obtained  by  the  convex- 
hull  construction.  Again,  we  have  computed  this  example  using  the  DG  method  with  a 
WENO  limiter  from  second  to  fourth  order  accuracy  and  have  obtained  similarly  good 
solutions,  even  though  only  the  second  order  result  is  shown  in  Figure  3.2.  The  high 
order  numerical  methods  are  observed  to  approximate  the  exact  solution  very  well  for 
this  example. 
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Figure  3.2:  Solid  lines:  the  exact  solution  of  (l.l)-(3.3)  with  the  initial  condition  (3.4)  at 
the  time  t  =  2;  solid  symbols:  the  numerical  solution  of  the  fifth  order  finite  volume  WENO 
scheme  (left)  and  the  second  order  discontinuous  Galerkin  method  with  WENO  limiters 
(right).  The  uniform  mesh  size  is  Ax  =  0.1. 

3.2  An  example  of  nonconvex  conservation  laws  with  poor  per¬ 
formance  of  high  order  schemes 

Despite  the  good  performance  of  high  order  numerical  methods  in  the  examples  demonstrated 
above,  there  are  also  examples  for  which  the  performance  of  high  order  numerical  methods 
is  poor.  In  this  subsection,  we  will  show  a  nonconvex  conservation  law,  such  that  many  high 
order  numerical  schemes  do  not  converge  to  the  unique  entropy  solution.  Indeed,  we  will 
show  in  the  following  that  the  numerical  solution  of  the  second  order  finite  volume  MUSCL 
scheme,  the  second  order  DG  method,  and  the  high  order  Godunov  type  scheme  with  fifth 
order  WENO  reconstruction  and  the  usual  time  step  governed  by  the  CFL  condition,  all 
with  the  Godunov  flux  (2.3),  would  stay  stationary  on  the  initial  data  for  a  specific  sequence 
of  meshes.  On  the  other  hand,  the  correct  entropy  solution  for  this  case  would  develop  a 
sequence  of  shocks  and  rarefaction  waves. 

Proposition  3.1.  Consider  the  semi-discrete  finite  volume  scheme  with  the  second  order 
MUSCL  reconstruction  (2.5)  and  the  monotone  Godunov  flux  (2.3),  for  solving  (1.1)  with 
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Figure  3.3:  Left:  the  nonconvex  flux  function  f(u)]  Right:  the  exact  solution  of  (l.l)-(3.5) 
with  the  initial  condition  (3.6)  at  the  time  t  =  2. 


the  nonconvex  flux  f(u)  defined  by 


fW) 


1, 

cos(57 r(y  —  1.8))  +  2.0, 
—  cos(57 r(y  —  2.2)), 

1, 


if  u  <  1.6 
if  1.6  <  u  <  2.0 
if  2.0  <  u  <  2.4 
if  u  >  2.4 


(see  the  left  panel  of  Figure  3.3)  and  the  initial  condition 


u0(x) 


1,  for  x  <  0 

3,  for  x  >  0 


(3.5) 


(3.6) 


The  numerical  solution  stays  stationary  over  the  time  evolution  for  a  specific  sequence  of 
meshes.  Therefore,  the  second  order  finite  volume  MUSCL  scheme  does  not  converge  to  the 
correct  entropy  solution  in  this  case  (the  right  panel  of  Figure  3.3). 


Proof:  In  order  to  prove  that  the  numerical  solution  stays  stationary,  it  is  sufficient  to  prove 
that  the  right  hand  side  of  (2.1)  equals  0  for  any  j. 

We  first  divide  the  computational  interval  into  N  =  (2 m  +  1)  evenly  spaced  subintervals 
such  that  the  origin  0  is  located  in  the  center  of  the  (m  +  l)-th  cell,  i.e.  xrn+\  =  0.  With  such 
a  specific  sequence  of  meshes,  and  under  the  initial  condition  (3.6),  the  initial  cell  averages 
are  given  by 


u°j  =  1,  for  j  =  1, ...,  m;  u°m+1  =  2; 


u°j  =  3,  for  j  =  m  +  2, ...,  N, 
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and  the  corresponding  reconstructed  point  values  by  the  MUSCL  procedure  (2.5)  at  the  cell 
boundaries  are 


=1,  for 

3+2 

J  = 

0, ...,  m  —  1; 

u±,0i  =3,  for 

3  +  2 

J  z 

0 

U  i  = 

m+2 

1, 

+,0  IK 

u  i  =  1.5, 

m+i 

u  0  3  =  2.5, 

m+i  ’ 

+,0 
u  i 

m+- 

If  we  use  the  monotone  Godunov  flux  (2.3)  as  the  numerical  flux,  then  fGi  =  1  for  all  j, 

2 

even  for  the  flux  at  xm+i  and  xmJ_3  where 


h m+i 


fG  i  =  fG(u  i ,u+  i )  =  /  T(l,  1.5)  =  max  f(u)  =  1, 

fG_  3  =  fG(u~  .  3,w+  ,  i)  =  fG( 2.5,3)  =  min  f(u)  =  1. 

Therefore,  the  right  hand  side  of  (2.1)  equals  0  for  all  j.  This  implies  that  the  cell  averages 
Uj ,  hence  the  numerical  solution  of  the  finite  volume  MUSCL  scheme,  stays  stationary  as 
(3.6)  over  time. 

The  correct  entropy  solution  can  be  solved  by  the  convex  hull  reconstruction.  The  entropy 
solution  at  the  time  t  —  2  is  plotted  in  the  right  panel  of  Figure  3.3.  The  discontinuity  at 
the  origin  in  the  initial  condition  generates  a  shock,  a  rarefaction  wave  followed  by  another 
shock  over  time.  Clearly  the  numerical  solution  of  the  finite  volume  MUSCL  scheme  fails  to 
converge  to  the  correct  entropy  solution  in  this  case.  Indeed,  we  can  easily  check  that  the 
stationary  numerical  solution  is  one  of  the  weak  solutions,  but  not  the  entropy  solution,  as 
it  violates  the  entropy  inequality  for  the  convex  entropies  U(u)  =  \u  —  c\  with  c  G  [2.0,  2.2], 
This  completes  the  proof.  ■ 


Proposition  3.2.  Consider  the  semi-discrete  DG  method  (2.11)  for  (1.1)  with  the  nonconvex 
flux  f(u)  defined  by  (3.5)  and  the  initial  condition  (3.6).  The  numerical  solution  of  the  second 
order  DG  method  with  the  Godunov  flux  (2.3)  will  stay  stationary  over  time  for  a  specific 
sequence  of  meshes.  Therefore,  the  DG  method  could  not  converge  to  the  correct  entropy 
solution  (the  right  panel  of  Figure  3.3). 
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Proof:  Let  the  sequence  of  meshes  be  the  same  as  that  in  the  proof  of  Proposition  3.1 


and  consider  the  second  order  DG  method  with  the  solution  space  and  test  space  to  be  the 
piecewise  linear  function  V^x  as  in  (2.8).  If  we  adopt  a  local  orthogonal  basis  over  Ij  as  {1, 
with  ,  the  semi-discrete  DG  method  is  formulated  as  finding  uAx  =  up  +  on 

Ij ,  such  that 


d  (o) 

mu> 

]_d_  (i) 

12  dtUj 


Axi 


fj+i  ~  fi- 1 

1 


f(uAx)dx  — 


Axj  \  Ax.j 


1  .  1  . 

9  fjH  +  9^-1 


Initially,  from  (3.6),  we  have 


1, 

for 

3 

=  1,2, 

"1 

m 

2, 

for 

J 

=  m  + 

1 

3, 

for 

3 

=  m  + 

2, 

...,1V 

o, 

for 

3 

=  1,2, 

*  *? 

m 

3, 

for 

3 

=  m  + 

1 

o, 

for 

3  = 

=  m  + 

2, 

...,N 

By  (2.3),  it  could  be  easily  checked  that  the  Godunov  flux 


(3.7) 

(3.8) 


(3.9) 


(3.10) 


ff  j  =  1- 


(3.11) 


and 


—  p  f(u  X)dx  =  1,  Vj, 


(3.12) 


even  for  the  flux  at  the  boundaries  of  the  (m  +  l)-th  cell.  Therefore,  for  any  j,  the  right  hand 
sides  of  both  (3.7)  and  (3.8)  vanish.  In  other  words,  both  the  mean  value  and  the  slope  of  the 
numerical  solution  over  each  cell  would  stay  stationary  over  time.  Therefore,  the  numerical 
solution  of  the  semi-discrete  DG  method  with  the  Godunov  flux  stays  stationary  with  respect 
to  the  initial  condition  (3.6)  and  could  not  converge  to  the  correct  entropy  solution  as  shown 
in  the  right  panel  of  Figure  3.3.  This  completes  the  proof.  ■ 


Proposition  3.3.  Consider  the  high  order  Godunov  type  scheme  (2.7)  with  the  fifth  order 
WENO  reconstruction  (2.6),  for  solving  (1.1)  with  the  nonconvex  flux  function  f{u)  defined 
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by  (3.5)  and  the  initial  condition  (3.6).  The  numerical  solution  stays  stationary  when  the 
time  step  is  small  enough  for  a  specific  sequence  of  meshes.  Numerical  solution  of  this  high 
order  Godunov  type  scheme  therefore  could  not  converge  to  the  correct  entropy  solution  (the 
right  panel  of  Figure  3.3). 

Proof:  The  proof  is  similar  to  that  in  [18],  and  is  thus  omitted.  ■ 

We  have  shown  that  the  numerical  solutions  for  (l.l)-(3.5)  with  the  initial  condition 
(3.6)  would  stay  stationary  for  many  high  order  schemes  including  the  second  order  finite 
volume  MUSCL  scheme,  the  second  order  discontinuous  Galerkin  method  and  the  fifth  order 
Godunov  type  WENO  scheme  with  the  Godunov  flux,  which  is  the  least  dissipative  among 
all  monotone  fluxes,  for  a  specific  sequence  of  meshes.  The  numerical  solutions  of  these 
schemes  with  more  dissipative  fluxes,  e.g.  the  Lax-Friedrichs  flux,  or  with  different  meshes, 
however,  may  not  necessarily  be  stationary.  In  the  following,  we  compare  the  performance 
of  different  schemes  with  different  types  of  numerical  fluxes  and  different  meshes  for  solving 
(1.1)  with  the  flux  f{u)  given  by  (3.5).  We  have  performed  extensive  numerical  experiments, 
but  are  presenting  only  a  selected  few  to  demonstrate  the  poor  performance  of  the  low  order 
monotone  schemes  (slow  convergence)  and  high  order  schemes  (slow  or  no  convergence). 


1.  Monotone  scheme  with  the  Lax-Friedrichs  flux  (the  left  panel  of  Figure  3.4)  and  with 
the  Godunov  flux  (the  right  panel  of  Figure  3.4)  for  the  initial  data  (3.6). 

2.  Monotone  scheme  with  the  Lax-Friedrichs  flux  (the  left  panel  of  Figure  3.5)  and  with 
the  Godunov  flux  (the  right  panel  of  Figure  3.5)  for  the  initial  condition 

3,  for  —  1  <  x  <  0 


Uq(x)  = 


1, 


for  0  <  x  <  1 


(3.13) 


and  a  periodic  boundary  condition. 


3.  Fifth  order  finite  volume  WENO  scheme  with  the  Lax-Friedrichs  flux  (the  left  panel 
of  Figure  3.6)  and  with  the  Godunov  flux  (the  right  panel  of  Figure  3.6)  for  the  initial 
data  (3.6). 
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4.  Discontinuous  Galerkin  method  with  WENO  limiters  (DG-WENO),  with  the  Lax- 
Friedrichs  flux  (the  left  panels  of  Figure  3.7)  and  with  the  Godunov  flux  (the  right 
panels  of  Figure  3.7)  for  the  initial  data  (3.6). 
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Figure  3.4:  The  numerical  solution  of  the  monotone  scheme  with  the  Lax-Friedrichs  flux 
(left)  and  the  Godunov  flux  (right)  for  the  nonconvex  scalar  conservation  law  (l.l)-(3.5) 
with  the  initial  condition  (3.6)  at  t  =  2  using  N  =  50  (unfilled  square)  and  N  =  100  (filled 
square)  uniform  cells. 


From  the  demonstrated  numerical  results,  we  can  make  the  following  observations. 

1.  The  first  order  monotone  scheme  with  the  Godunov  flux  has  a  much  better  performance 
than  the  one  with  the  Lax-Friedrichs  flux  (see  Figure  3.4).  However,  with  a  periodic 
boundary  condition,  after  shocks  and  rarefaction  waves  interact  with  each  other,  the 
convergence  is  very  slow  even  for  the  first  order  Godunov  scheme  (see  Figure  3.5). 

2.  The  fifth  order  finite  volume  WENO  scheme  with  the  Godunov  flux  might  not  converge 
to  the  correct  entropy  solution.  The  scheme  with  the  Lax-Friedrichs  flux  seems  to 
converge  to  the  entropy  solution  with  a  slow  convergence  rate  (the  left  panel  of  Figure 
3.6),  which  might  be  related  to  the  fact  that  the  reconstruction  of  the  solution  at  the 
rarefaction  wave  comes  from  neighboring  cells  and  is  not  a  good  approximation  when 
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Figure  3.5:  The  numerical  solution  of  the  monotone  scheme  with  the  Lax-Friedrichs  flux 
(left)  and  the  Godunov  flux  (right)  for  the  nonconvex  scalar  conservation  law  (l.l)-(3.5) 
with  the  initial  condition  (3.13)  ait  —  2  with  a  periodic  boundary  condition  using  N  =  200 
(unfilled  square),  N  =  1000  (filled  square)  and  N  =  36000  (solid  line)  uniform  cells. 

the  rarefaction  wave  is  surrounded  by  two  shocks  at  its  early  stage  of  development. 
Further  numerical  experiments  on  (l.l)-(3.5)  with  the  initial  condition  (3.13)  indicate 
that  the  finite  volume  WENO  scheme  with  the  Lax-Friedrichs  flux  could  not  converge 
to  the  entropy  solution. 

3.  The  DG  scheme  with  WENO  limiters  might  not  converge  to  the  correct  entropy  solution 
of  (1.1).  Different  meshes,  fluxes  and  solution  spaces  in  the  DG  scheme  will  give 
different  numerical  results.  The  numerical  solutions  seem  to  converge  to  the  correct 
entropy  solution  in  some  but  not  all  cases. 

3.3  Another  nonconvex  conservation  law:  f{u)  =  sin(u) 

In  section  3.2,  we  have  shown  a  nonconvex  conservation  law,  with  which  the  numerical  so¬ 
lution  of  the  high  order  schemes  could  not  converge  to  the  correct  entropy  solution.  Indeed 
this  is  not  the  only  example.  We  will  show  in  the  following  another  example  of  poor  perfor¬ 
mance  of  high  order  numerical  methods.  Consider  the  Riemann  problem  of  the  nonconvex 
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Figure  3.6:  The  numerical  solution  of  the  fifth  order  finite  volume  WENO  scheme  with  the 
Lax-Friedrichs  flux  (left)  and  the  Godunov  flux  (right)  for  the  nonconvex  scalar  conservation 
law  (l.l)-(3.5)  with  the  initial  condition  (3.6)  at  t  —  2  using  N  —  50  (unfilled  square)  and 
N  =  800  (filled  square)  uniform  cells. 

conservation  law  (1.1)  with  the  flux  function 


f(u)  =  sin(w), 


and  the  initial  condition 


f  7t/64,  if  x  <  0, 

\  255tt/64,  if  x  >  0. 


(3.14) 


It  is  shown,  in  Figure  3.8  for  the  high  order  finite  volume  schemes  and  in  Figure  3.9  for  the 
discontinuous  Galerkin  methods,  that  the  numerical  solution  of  these  high  order  schemes 
would  not  always  converge  to  the  entropy  solution.  One  of  the  rarefaction  waves  in  the 
compound  wave  is  sometimes  missing  in  the  solutions  of  these  high  order  schemes. 
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Figure  3.7:  The  numerical  solution  of  the  DG-WENO  schemes  at  t  —  2  for  (l.l)-(3.5)  with 
the  initial  condition  (3.6)  and  the  Lax-Friedrichs  flux  (left)  and  the  Godunov  flux  (right) 
using  N  —  51  (unfilled  square)  and  N  =  101  (filled  square)  uniform  cells.  From  top  to 
bottom  are  figures  for  the  second,  third  and  fourth  order  schemes  respectively. 


X  X 

Figure  3.8:  Solid  lines:  the  exact  solution  of  (1.1)  with  f(u)  =  sin(-u)  and  the  initial  condition 
(3.14)  at  the  time  t  =  4;  solid  symbols:  the  numerical  solution  of  the  second  order  MUSCL 
scheme  (left),  and  the  fifth  order  finite  volume  WENO  scheme  (right).  The  numerical  flux 
is  the  Godunov  flux.  The  mesh  is  uniform  with  Ax  =  0.05. 


4  First  order  monotone  modification 


In  this  section,  we  propose  a  first  order  monotone  modification  for  high  order  numerical 
methods.  The  objective  is  to  maintain  high  order  accuracy  in  smooth  regions  and  enforce 
convergence  to  the  entropy  solution  for  general  nonconvex  conservation  laws.  The  main¬ 
tenance  of  high  order  accuracy  is  achieved  by  a  carefully  designed  discontinuity  indicator 
[16]. 


4.1  Discontinuity  indicator 


In  order  to  maintain  the  high  order  accuracy  in  smooth  regions,  we  use  the  discontinuity 
indicator  designed  in  [16].  Specifically,  the  discontinuity  indicator  cpj  is  defined  as 


where 


(4.1) 


=|  Uj_ i  -  Uj  \2  +£,  =|  Uj+ 1  -  Uj- 1  | J  +£,  Pj  = 


Zj 


+ 

OLj- 1  CXj+ 2 


7 j  = 


( Ur 


—  Ur 


a.j 
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Figure  3.9:  Solid  lines:  the  exact  solution  of  (1.1)  with  f(u)  =  sin(w)  and  the  initial  condition 
(3.14)  at  the  time  t  —  4;  solid  symbols:  the  discontinuous  Galerkin  method  with  WENO 
limiters,  of  second  order  accuracy  (left),  third  order  accuracy  (middle)  and  fourth  order 
accuracy  (right).  The  numerical  flux  is  the  Godunov  flux.  The  mesh  is  uniform  with  Ax  = 
0.05. 


Here  Uj  refers  to  the  cell  average  of  the  numerical  solution  on  Ij  in  the  finite  volume  scheme 
and  the  discontinuous  Galerkin  method,  £  is  a  small  positive  number  taken  as  10~6  in  the 
code,  and  umax  and  ur nm  are  the  maximum  and  minimum  values  of  Uj  over  all  cells.  The 
discontinuity  indicator  <j)j  has  the  property  that 

•  0  <  fa  <  1. 

•  (f)j  is  on  the  order  of  0(  Ax2)  in  smooth  regions. 

•  fa  is  close  to  0(1)  near  a  strong  discontinuity. 

4.2  Modification  on  the  fifth  order  finite  volume  WENO  scheme 


In  this  section,  we  propose  a  first  order  modification  to  the  fifth  order  finite  volume  WENO 
scheme.  Recall  from  section  2.1  that  the  cell  average  Uj  in  the  finite  volume  WENO  scheme 
is  updated  by 


duj 

dt 


(/i+* 


Axj  v  j+* 


(4.2) 


m 


which  the  numerical  flux  fJ+i  =  f(Uj+  or  more  precisely,  the  reconstructed  point 


values  Uj+i  at  the  cell  boundaries  play  an  essential  role.  The  modification  proposed  is  based 
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on  modifying  u±  k  at  the  cell  boundaries.  The  scheme  can  be  summarized  as  following,  after 
■?+2 

a  suitable  initialization  to  obtain  u°. 


1.  Perform  the  WENO  reconstruction. 

At  each  cell  interface,  say  x-+i,  reconstruct  the  point  values  vfl  i  using  neighboring 

2  J  i  2 

cell  averages  un  by  the  fifth  order  WENO  reconstruction  procedure  described  in  section 
2.1. 

2.  Identify  the  troubled  cell  boundary  Xj+ 1. 

Criterion  I:  A  cell  boundary  is  good,  if  u^+1,  Uj  and  iij+ 1  all  fall  into  the  same 
linear,  convex  or  concave  region  of  the  flux  function  f(u).  Otherwise,  it  is  defined  to 
be  a  troubled  cell  boundary. 


When  the  flux  function  f[u)  is  nonconvex,  one  can  divide  its  domain  into  several 
subregions,  across  each  of  which  f[u)  is  purely  linear,  convex  or  concave.  For  example, 
for  f[u)  defined  in  (3.5),  there  are  several  linear  regions,  e.g.  [— oo,  1.6],  [2.4,  oo],  convex 
regions,  e.g.  [1.6,  1.7],  [1.9,  2.0],  [2.1,  2.3],  and  concave  regions,  e.g.  [1.7,  1.9],  [2.0, 
2.1],  [2.3,  2.4].  Points  are  said  to  fall  into  the  same  region  if  all  of  them  are  contained 
in  the  same  subregion. 


3.  At  troubled  cell  boundaries,  modify  the  numerical  flux  fJ+i  with  a  discontinuity  indi¬ 
cator. 


Let  fj+ 1  =  where 


W+i  =  (!  -  ,  I  +  u =  (!  -  +  <f>jUj+ i, 


j+i 


3'  j+ 


j+i 


(4.3) 


with  <fj  defined  by  (4.1),  if  Xj+ 1  is  a  troubled  cell  boundary.  Otherwise,  at  good  cell 

boundaries,  um;f  =  u±,  1 . 

J+2  ^+2 

4.  Evolve  the  cell  averages  u3  by  (4.2). 
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Remark.  When  a  troubled  cell  boundary  is  at  a  strong  discontinuity,  4>j  ~  1,  hence  ~ 
Uj  and  ~  Uj+ 1,  indicating  a  first  order  monotone  scheme  is  taking  the  effect  at  a 

nonconvex  discontinuity  region.  When  a  troubled  cell  boundary  is  in  a  smooth  region,  the 
modification  is  obtained  with  the  magnitude  at  most  of  the  size 


max 


M, 


U 


uj+ 1 


u 


0{  Arc5), 


hence  it  does  not  affect  the  fifth  order  accuracy  of  the  scheme. 

Numerical  examples  are  given  below  to  check  the  maintained  order  of  accuracy  for  smooth 
solutions  and  convergence  towards  the  entropy  solution. 


Example  4.1.  The  nonconvex  conservation  law 


u " 


Ut  +  (  y  )  =0,  u0(x)  =  sin(Tnr). 


(4.4) 


Table  4.1  gives  the  L 1  errors  and  the  corresponding  orders  of  accuracy  of  the  regular  and 
modified  finite  volume  WENO  scheme.  The  L 1  errors  of  both  schemes  are  comparable.  Very 
little  difference  is  observed. 


Table  4.1:  The  L 1  errors  and  the  corresponding  orders  of  accuracy  for  the  regular  and 
modified  finite  volume  WENO  schemes  with  the  Godunov  flux  for  (4.4)  at  the  time  t  =  0.2. 


No.  of 

regular 

modified 

points 

L 1  error 

order 

L 1  error 

order 

100 

2.128E-05 

2.128E-05 

200 

1.120E-06 

4.25 

1.120E-06 

4.25 

300 

1.648E-07 

4.73 

1.648E-07 

4.73 

400 

4.079E-08 

4.85 

4.079E-08 

4.85 

500 

1.366E-08 

4.90 

1.366E-08 

4.90 

Example  4.2.  The  nonconvex  conservation  law  (1.1)  with  f(u)  =  sin(u)  and  the  initial 
condition  (3.14).  The  solutions  at  the  time  t  =  4  are  plotted  in  the  left  panel  of  Figure 
4.1.  As  shown,  the  numerical  solutions  of  the  modified  scheme  successfully  converge  to  the 
correct  entropy  solution,  with  the  development  of  the  complex  solution  structure  containing 
a  shock,  a  rarefaction  wave,  followed  by  another  shock  and  another  rarefaction  wave. 
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Figure  4.1:  Left  figure:  Example  4.2.  Solid  lines:  the  exact  solution  at  time  t  —  4;  solid 
symbols:  the  numerical  solution  of  the  modified  fifth  order  finite  volume  WENO  scheme  with 
the  uniform  mesh  size  Ax  =  0.1  (unfilled  squares),  Ax  =  0.05  (solid  symbols).  Right  figure: 
Example  4.3.  Solid  lines:  the  reference  solution  (computed  by  the  first  order  monotone 
Godunov  scheme  with  N  =  36000  points)  at  time  t  =  2;  symbols:  the  numerical  solution 
of  the  modified  fifth  order  finite  volume  WENO  scheme  with  N  =  200  (unfilled  square)  and 
N  =  1000  (filled  square)  uniform  cells.  The  numerical  flux  in  both  figures  is  the  Godunov 
flux. 

Example  4.3.  The  nonconvex  conservation  law  (1.1)  with  f(u)  given  by  (3.5),  the 
initial  condition  (3.13)  and  a  periodic  boundary  condition.  As  shown  in  the  right  panel  of 
Figure  4.1,  the  numerical  solution  of  the  modified  scheme  converges  to  the  correct  entropy 
solution,  with  much  faster  convergence  rate  when  compared  to  the  first  order  monotone 
scheme  (Figure  3.5).  We  also  test  our  modified  scheme  by  further  refining  the  numerical 
meshes.  Convergence  to  the  entropy  solution  is  observed. 

4.3  Modification  on  the  discontinuous  Galerkin  methods 

In  this  subsection,  we  propose  a  first  order  modification  to  the  high  order  discontinuous 
Galerkin  methods.  Recalling  from  section  2.2,  the  piecewise  polynomial  over  each  cell  Ij  is 
updated  by  its  moments  uj\t)  as  defined  in  (2.10)  through  (2.11).  To  modify  high  order 
DG  methods,  we  make  the  following  two  observations: 

1.  Only  the  numerical  flux  fJ+i ,  or  more  precisely,  the  left  and  right  limits  of  discontinuous 
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solution  u  i  and  u++i  at  the  cell  interface  are  involved  in  updating  the  cell  average 

(o) 

u)>. 

2.  The  whole  polynomial  over  Ij  is  involved  in  updating  the  higher  order  moments 
for  l  >  1  through  the  integration  term  on  the  right  hand  side  of  (2.11). 


Based  on  these  two  observations,  a  monotone  modification  for  (k  +  l)-th  (k  >  1)  order 
DG  methods  with  WENO  limiters  is  designed  as  follows.  After  a  suitable  initialization, 

uf*  —  —  [  u(x,0)vl'j\x)dx,  l  =  0, 1,  2, ...,  k. 

('i  ■' h 

1.  Apply  the  WENO  limiters  in  the  oscillatory  troubled  cells. 

The  oscillatory  troubled  cells  refer  to  those  cells  where  oscillations  at  discontinuities 
might  occur  in  high  order  DG  methods.  We  refer  to  [12]  for  the  details  of  identifying 
such  oscillatory  troubled  cells  and  reconstructing  high  order  polynomials  there  by  the 
WENO  procedure. 

2.  Identify  the  nonconvex  troubled  cells. 

Criterion  II:  A  cell  Ij  is  called  a  good  cell,  if  the  numerical  solution  uAx  in  the  cells 
Ij- 1  ,  Ij  and  Ij+ 1  falls  into  the  same  linear,  convex  or  concave  region  of  the  flux  function 
f(u).  Otherwise,  it  is  defined  to  be  a  nonconvex  troubled  cell. 

3.  Evolve  the  polynomial  solution. 


If  a  cell  Ij  is  good,  then  evolve  the  polynomial  function  vAx | j.  by  the  regular  DG 
method  through  (2.11). 


Otherwise,  evolve  the  cell  average  ra  f}  only  by 


du 


(o) 


dt 


A  Xj  i+5 


where  the  numerical  flux  /-,  i  =  fimfi ,  umf\)  with 

J  2  J  i  9  J~r  9 


(4.5) 


A2..(0) 
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with  (f)j  defined  by  (4.1)  at  the  troubled  cell  boundaries  as  defined  via  Criteria  I 
in  section  4.2,  and  um’ f  =  1  at  good  cell  boundaries. 

4.  Reconstruct  the  polynomials  in  the  troubled  cells. 

For  a  troubled  cell  Ij,  where  only  the  cell  average  is  updated,  the  polynomial  is  recon¬ 
structed  from  neighboring  cell  averages  in  a  WENO  fashion,  as  in  [12], 


Example  4.4.  The  nonconvex  conservation  law  (4.4).  Tables  4.2  and  4.3  give  the  L 1 
errors  and  the  corresponding  orders  of  accuracy  of  the  regular  DG-WENO  methods  and  the 
modified  DG-WENO  methods  respectively.  The  tables  show  comparable  L 1  errors. 


Table  4.2:  The  L 1  errors  and  the  corresponding  orders  of  accuracy  for  the  regular  DG-WENO 
method  with  the  Godunov  flux  for  (4.4)  at  the  time  t  =  0.2. 


No.  of 

k  =  1 

k  =  2 

k  =  3 

points 

L1  error 

order 

L1  error 

order 

Ll  error 

order 

30 

5.207E-03 

7.655E-04 

3.679E-04 

60 

1.230E-03 

2.08 

4.889E-05 

3.97 

5.403E-06 

6.09 

90 

5.365E-04 

2.05 

1.265E-05 

3.33 

6.844E-07 

5.10 

120 

2.862E-04 

2.18 

5.122E-06 

3.14 

2.045E-07 

4.20 

150 

1.648E-04 

2.47 

2.608E-06 

3.02 

8.341E-08 

4.02 

Table  4.3:  The  L 1  errors  and  the  corresponding  orders  of  accuracy  for  the  modified  DG- 
WENO  method  with  the  Godunov  flux  for  (4.4)  at  the  time  t  =  0.2. 


No.  of 

k  =  1 

k  =  2 

k  =  3 

points 

L 1  error 

order 

L 1  error 

order 

L 1  error 

order 

30 

2.313E-03 

1.146E-03 

6.951E-04 

60 

1.242E-03 

4.40 

4.924E-05 

4.54 

5.714E-06 

6.93 

90 

5.382E-04 

2.06 

1.264E-05 

3.35 

7.044E-07 

5.16 

120 

2.866E-04 

2.19 

5.120E-06 

3.14 

2.048E-07 

4.29 

150 

1.649E-04 

2.48 

2.608E-06 

3.02 

8.291E-08 

4.05 

Example  4.5.  The  nonconvex  conservation  law  (1.1)  with  f(u)  =  sin(u)  and  the  initial 
condition  (3.14).  Its  solution  at  the  time  t  =  4  of  the  second  and  fourth  order  DG-WENO 
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methods  are  plotted  in  Figure  4.2.  As  shown,  the  numerical  solutions  of  the  modified  schemes 
successfully  converge  to  the  correct  entropy  solution. 


X  X 

Figure  4.2:  Solid  lines:  the  exact  solution  of  (1.1)  with  (3.14)  at  the  time  t  —  4;  solid 
symbols:  the  modified  discontinuous  Galerkin  method  with  WENO  limiters,  of  second  order 
accuracy  (left)  and  fourth  order  accuracy  (right).  The  numerical  flux  is  the  Godunov  flux. 
The  uniform  mesh  sizes  are  Ax  =  0.1  (unfilled  squares)  and  Ax  =  0.05  (solid  symbols), 
respectively. 


Example  4.6.  The  nonconvex  conservation  law  (1.1)  with  f{u)  defined  in  (3.5),  the 
initial  condition  (3.13)  and  a  periodic  boundary  condition.  As  shown  in  Figure  4.3,  the 
numerical  solutions  of  the  modified  schemes  converge  to  the  correct  entropy  solution  with 
much  faster  convergence  rate  when  compared  with  the  first  order  monotone  scheme.  We 
also  test  our  modified  scheme  by  further  refining  the  numerical  mesh.  Convergence  to  the 
entropy  solution  is  observed. 

5  Second  order  modification  with  an  entropic  projec¬ 
tion 

A  MUSCL  type  method  with  an  entropic  projection,  enjoying  the  cell  entropy  inequality 
for  all  convex  entropy  functions,  is  proposed  in  [1],  This  entropy  satisfying  property  for  all 
convex  entropy  functions  is  very  attractive  in  the  design  of  high  order  numerical  methods 
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Figure  4.3:  Solid  lines:  the  reference  solution  obtained  by  the  first  order  monotone  Godunov 
scheme  with  N  =  36000  grid  points  for  (1.1)  with  (3.5),  the  initial  condition  (3.13)  and 
a  periodic  boundary  condition  at  the  time  t  =  2.  Symbols:  the  modified  discontinuous 
Galerkin  method  with  WENO  limiters,  of  second  order  accuracy  (left),  third  order  accuracy 
(middle)  and  fourth  order  accuracy  (right)  with  N  =  200  (unfilled  squares)  and  N  =  1000 
(solid  squares)  uniform  cells.  The  numerical  flux  is  the  Godunov  flux. 


for  general  nonconvex  conservation  laws.  Based  on  the  idea  of  this  entropic  projection,  we 
propose  a  second  order  modification  for  the  high  order  schemes. 

5.1  Review  of  the  MUSCL  method  satisfying  all  the  numerical 
entropy  inequalities 

The  scheme  in  [1]  uses  piecewise  linear  functions  as  the  solution  space.  Specifically,  the 
numerical  solution  at  time  level  n  can  be  written  as  un  =  u™  +  with  ^  =  ‘^-L  over  the 
cell  Ij.  It  consists  of  two  steps  to  evolve  from  un  to  un+1. 

1.  Exact  evolution  (T&t):  Evolve  (1.1)  exactly  for  a  time  step  At,  to  obtain  a  solution 
un+1,  which  in  general  is  not  a  piecewise  linear  function  anymore. 

2.  An  entropic  projection  ( P x):  Find  a  second  order  approximation  to  un+l  by  a  piecewise 
linear  function  un+1 ,  satisfying 

[  U(un+1(x))dx  <  [  U(un+1(x))dx,  Vj 

Jij  Jij 
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(5.1) 


for  all  convex  entropy  functions  U(u).  Second  order  reconstruction  satisfying  (5.1)  can 
be  obtained  by  setting  the  cell  average  as 


un+l  = _ 

i  Ax 


u 


n+ 1 


j  Jh 


(5.2) 


and  the  slope  as 


where 


s™+1  =  Dun+1  |/.=  minmod  i((y) 


(5.3) 


<(v)  = 


rx  i 

^  un+l(x)dx  — 


un+1(x)dx 


y~xi-h 


AxJ  \^+±  -yjy 
The  minmod  function  of  g{x )  on  the  interval  (a,  b)  is  defined  as 


0, 


if  32/1,2/2  e  (a,  5),  s.t.  g(yi)g(y2)  <  0, 


mmmodM9(z)  =  {  SS9'9’’  *  9<9>  >  °-  V«  £  <«-»>- 

max  (7(1/),  if  ^(y)  <0,  Vy  G  (a,  6). 

(a, 6) 

In  summary,  the  scheme  can  be  written  out  in  the  following  abstract  form 


(5.4) 


(5.5) 


«,*«  =  pi  0  TM(un)  =  Q1(At)(un). 


(5.6) 


It  enjoys  the  following  convergence  theorem  as  proved  in  [1],  Comparing  with  the  high 
order  Godunov  type  schemes  (2.7)  which  may  fail  to  satisfy  at  least  some  entropy  conditions 
[18],  the  entropic  projection  (P1),  especially,  the  reconstruction  of  the  slope  s"+1  in  (5.3),  is 
crucial  to  ensure  convergence  towards  the  entropy  solution. 


Theorem  5.1  [1].  Let  T  =  nAt,  u(-,T)  be  the  exact  entropy  solution  to  (1.1)  with  the 
initial  data.  u0  G  L 1  r)BV(lZ),  /qq  =  maxue[minU0jmaxU0]  f(u),  then  there  exists  a  constant  C, 
such  that 

\\Q\At)nu°  -u(-,T)\\Lr  <C(fooVT^t  +  Ax^T/At).  (5.7) 

Therefore,  the  second  order  MUSCL  scheme  with  the  entropic  projection  (5.6)  converges  to 
the  unique  entropy  solution. 
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5.2  Second  order  schemes  with  the  entropic  projection 

In  this  subsection,  the  entropic  projection  (P1),  or  more  explicitly,  the  reconstruction  of 
the  slope  s”+1  in  (5.3),  is  used  to  design  a  second  order  method  for  general  nonconvex 
conservation  laws.  With  a  suitable  initialization  to  obtain  a  piecewise  linear  function,  u®  + 
s^j  with  the  scheme  consists  of  the  following  steps  in  each  time  step  evolution. 

1.  Identify  nonconvex  troubled  cells,  for  which  we  refer  to  Criterion  II  in  section  4.3  for 
the  details. 

2.  Update  u™+1  for  the  good  cells. 

If  a  cell  Ij  is  good,  evolve  u™  by  the  regular  finite  volume  method  through  (4.2), 
where  f\,i  =  f(uZ  i ,  u+  C  with  uZ  i  and  u+  1  being  the  left  and  right  limits  of  the 
numerical  solution  at  the  cell  interface. 

3.  Update  uj+1  and  s"+1  in  nonconvex  troubled  cells. 

If  a  cell  Ij  is  a  nonconvex  troubled  cell,  update  un+1\ij  by  (5.2)  and  (5.3),  where 
un+1\ij  is  approximated  by  a  first  order  refined  mesh  evolution.  Specifically,  a  first 
order  monotone  scheme  is  used  to  numerically  evolve  (1.1)  with  the  initial  condition 

Ui+Si£i,  for  x  e  Ii,  l=j-l,j,j+l  (5-8) 

and  a  periodic  boundary  condition  on  Ij_i  U  Ij  U  IJ+\  for  time  At.  We  remark  that 
periodic  boundary  condition  is  allowed  here  as  the  information  outside  the  domain, 
Ij- 1  U  Ij  U  Ij+ 1,  would  not  affect  the  solution  on  Ij  after  time  At,  which  is  restricted 
by  the  CFL  condition.  Let  Ij  be  uniformly  discretized  by 

Ij  m.=l[y  m— ^  i  it  m.+  i],  Vm+ ^  Vm—  i  Ax/N , 

then  «n+1|j  is  approximated  by  a  piecewise  constant  function  sitting  on  the  refined 
numerical  mesh  with  the  truncation  error  ~  O(Sx)  =  0( Ax2).  (5.3)  is  numerically 
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implemented  by  finding  y  among  the  (finitely  many)  refined  cell  boundaries  ym+i 
(rn  =  0, N )  that  minimizes  (( y )  as  defined  in  (5.4). 

4.  Update  s"+1  in  the  good  cells. 

In  the  good  cells,  reconstruct  the  slope  s"+1  in  the  MUSCL  fashion  from  neighboring 
cell  averages. 


Remark.  The  implementation  of  the  entropic  projection  proposed  above  is  very  compu¬ 
tational  expensive  and  impractical.  An  efficient  implementation  for  specific  equations,  for 
example  for  the  Burgers  equation  with  the  flux  function  f(u)  =  ^  is  available  in  [1],  For 
general  nonconvex  flux  functions,  an  efficient  implementation  for  the  entropic  projection  is 
nontrivial  and  is  worth  a  further  investigation. 

In  the  following,  we  provide  numerical  examples  of  the  proposed  second  order  scheme 
with  the  entropic  projection. 

Example  5.1.  The  nonconvex  conservation  law  (4.4).  Table  5.1  gives  the  L 1  errors  and 
the  corresponding  orders  of  accuracy  for  the  regular  second  order  finite  volume  scheme  with 
the  MUSCL  reconstruction  and  the  proposed  scheme  with  the  entropic  projection.  The  table 
shows  comparable  L 1  errors  and  second  order  accuracy  as  expected. 


Table  5.1:  The  L 1  errors  and  the  corresponding  orders  of  accuracy  for  the  regular  second 
order  finite  volume  scheme  with  the  MUSCL  reconstruction  and  the  proposed  scheme  with 
the  entropic  projection  for  (4.4)  at  the  time  t  =  0.2. 


No.  of 

MUSCL 

entropic  projection 

points 

L1  error 

order 

L 1  error 

order 

60 

3.198E-03 

3.202E-03 

120 

8.863E-04 

1.85 

8.865E-04 

1.85 

180 

4.129E-04 

1.88 

4.129E-04 

1.88 

240 

2.394E-04 

1.90 

2.394E-04 

1.90 

300 

1.548E-04 

1.95 

1.548E-04 

1.95 

360 

1.095E-04 

1.90 

1.095E-04 

1.90 
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Example  5.2.  The  nonconvex  conservation  laws  (1.1)  with  f(u)  =  sin(w)  and  the  initial 
condition  (3.14).  The  solutions  at  the  time  t  —  4  are  plotted  in  the  left  panel  of  Figure  5.1. 
As  shown,  the  numerical  solution  of  the  scheme  with  the  entropic  projection  successfully 
converges  to  the  correct  entropy  solution,  with  better  resolution  of  the  rarefaction  wave 
when  compared  with  the  regular  MUSCL  reconstruction  (the  left  panel  of  Figure  3.8). 
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Figure  5.1:  Left  figure:  the  solid  line  is  the  exact  solution  of  (1.1)  with  (3.14)  at  the  time 
t  —  4,  while  the  solid  symbols  are  the  numerical  solution  of  the  second  order  finite  volume 
scheme  with  the  entropic  projection  with  Ax  =  0.1  (unfilled  squares),  Ax  =  0.05  (solid 
symbols).  Right  figure:  the  solid  line  is  the  reference  solution  (obtained  with  a  first  order 
monotone  Godunov  scheme  with  N  =  36000  grid  points)  at  the  time  t  —  2,  while  the 
symbols  are  the  numerical  solution  of  the  second  order  scheme  with  the  entropic  projection 
with  N  =  200  (unfilled  squares)  and  N  =  1000  (solid  symbols)  uniform  cells.  The  numerical 
flux  in  both  figures  is  the  Godunov  flux. 


Example  5.3.  The  nonconvex  conservation  law  (1.1)  with  f(u)  defined  by  (3.5),  the 
initial  condition  (3.13)  and  a  periodic  boundary  condition.  As  shown  in  the  right  panel  of 
Figure  5.1,  the  numerical  solution  of  the  second  order  scheme  with  the  entropic  projection 
converges  to  the  correct  entropy  solution,  with  much  faster  convergence  rate  when  compared 
with  the  first  order  monotone  scheme  (Figure  3.5). 
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5.3  Second  order  modification  to  the  fifth  order  finite  volume 
WENO  schemes 


In  this  subsection,  the  second  order  method  just  designed  is  used  as  a  building  block  to 
modify  the  fifth  order  finite  volume  WENO  scheme.  The  modification  follows  a  similar  line 
as  the  first  order  monotone  modification  described  in  section  4.2.  The  procedure  is  outlined 
as  following. 


At  each  time  step  evolution,  we  would  like  to  update  Uj  =  {uv  up  ?i)}  over  the  cell 


Ij.  At  the  initial  stage,  is  obtained  by  the  WENO  reconstruction  from  u.  ul  and  urj 
refer  to  approximations  to  the  left  and  right  boundaries  of  L.  u\  =  u+_  1  and  ur.  =  u" 


0  2 


3  i+i 


1.  Identify  the  troubled  cell  boundaries,  for  which  we  refer  to  Criterion  I  in  section  4.2 
for  the  details. 


2. 


Modify  the  numerical  flux  fJ+i 
Specifically,  let  fj+i  =  f(u™f  7, 


with  a  discontinuity  indicator. 
u'n’t ),  where 


+  <Lt  =  (1  -  tfXi 


i+i 


j+2 


+  4>)u 


i  i+i’ 


(5.9) 


with  4>j  defined  by  (4.1)  at  the  troubled  cell  boundary. 


m,d=  ± 

U3+k  -  U3+¥ 


At  good  cell  boundaries, 


3.  Identify  nonconvex  troubled  cells. 

Criterion  IV:  A  cell  Ij  is  called  a  good  cell,  if  um  =  {um,  ulm,  u^}  with  m  = 
j  —  1  ,j,j  +  1,  fall  into  the  same  linear,  convex  or  concave  region  of  the  flux  function 
f(u).  Otherwise,  it  is  defined  to  be  a  nonconvex  troubled  cell. 


4.  Update  ul3  and  urj  for  a  nonconvex  troubled  cell  Ij. 

Perform  a  first  order  refined  mesh  evolution.  Specifically,  the  first  order  monotone 
scheme  is  used  to  numerically  evolve  (1.1)  with  the  initial  condition  (5.8)  where 


sm  =  2minmod(u(n  -  um,  um  -  ulm ),  m  =  j-  l,j,j  +  l 
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and  a  periodic  boundary  condition  to  obtain  un+1\jp  the  details  of  this  step  is  given  in 
section  5.2.  un+1  is  then  approximated  by  un-+1  +  s™+1£j  with  un-+l  obtained  from  (5.2) 
and  s”+1  from  (5.3). 


,n+l  _  ~n+ 1  _  \  n+1 
3  U3  2  j  ’ 


r,n+ 1 
U- 


=  u]+1  + 


1  „n+i 

2 


5.  Update  the  cell  averages  un+1  by  (4.2). 

6.  Update  i  by  the  WENO  reconstruction  procedure. 

7.  For  a  good  cell  Ij,  update  ul3  and  urj  by  setting  ulj  =  u+_1  and  urj  =  u~+1 


Remark.  Comparing  with  the  monotone  modification  of  high  order  schemes,  there  are  two 
additional  pieces  of  information  in  each  time  step  evolution,  namely  ulj  and  w).  They  are 
second  order  approximations  to  the  left  and  right  boundaries  of  a  troubled  cell  Ij  and  are 
used  in  (5.9)  to  modify  the  numerical  flux. 

In  the  following,  we  provide  numerical  examples  of  the  modified  fifth  order  finite  volume 
scheme  with  the  second  order  entropic  projection. 


Example  5.4.  The  nonconvex  conservation  law  (4.4).  Table  5.2  shows  comparable  L 1 
errors  and  the  corresponding  orders  of  accuracy  for  the  regular  fifth  order  finite  volume 
scheme  and  the  modihed  one  with  the  second  order  entropic  projection. 


Table  5.2:  The  L 1  errors  and  the  corresponding  orders  of  accuracy  for  the  regular  fifth  order 
finite  volume  WENO  scheme  and  the  modihed  one  with  the  second  order  entropic  projection 
for  (4.4)  at  the  time  t  =  0.2.  The  numerical  flux  is  the  Godunov  flux. 


No.  of 

regular 

modihed 

points 

L1  error 

order 

L1  error 

order 

100 

2.128E-05 

2.128E-05 

200 

1.120E-06 

4.25 

1.120E-06 

4.25 

300 

1.648E-07 

4.73 

1.648E-07 

4.73 

400 

4.079E-08 

4.85 

4.079E-08 

4.85 

500 

1.366E-08 

4.90 

1.366E-08 

4.90 

34 


Example  5.5.  The  nonconvex  conservation  law  (1.1)  with  f(u)  =  sin(ii)  and  the  initial 
condition  (3.14).  The  solutions  at  time  t  —  4  are  plotted  in  the  left  panel  of  Figure  5.2. 
As  shown,  the  numerical  solution  of  the  modified  fifth  order  finite  volume  scheme  with  the 
entropic  projection  successfully  converges  to  the  correct  entropy  solution. 


Figure  5.2:  Left:  Example  5.5.  The  solid  line  is  the  exact  solution  of  (1.1)  with  (3.14)  at  the 
time  t  —  4,  while  the  solid  symbols  are  the  numerical  solution  of  the  modified  fifth  order  finite 
volume  scheme  with  the  entropic  projection.  The  mesh  size  Ax  =  0.1  (unfilled  squares)  and 
Ax  =  0.05  (solid  symbols).  Right:  Example  5.6.  Solid  lines:  the  reference  solution  obtained 
from  the  first  order  monotone  Godunov  scheme  with  N  =  36000  grid  points  at  the  time 
t  —  2;  symbols:  the  numerical  solution  of  the  modified  fifth  order  finite  volume  WENO 
scheme  with  the  entropic  projection  with  N  =  200  (unfilled  square)  and  N  =  400  (filled 
square)  uniform  cells.  The  numerical  flux  in  both  figures  is  the  Godunov  flux. 


Example  5.6.  The  nonconvex  conservation  law  (1.1)  with  f(u)  defined  by  (3.5),  the 
initial  condition  (3.13)  and  a  periodic  boundary  condition.  As  shown  in  the  right  panel  of 
Figure  5.2,  the  numerical  solution  of  the  fifth  order  finite  volume  WENO  scheme  with  the 
second  order  entropic  projection  approximates  the  correct  entropy  solution  very  well  with 
relatively  coarse  meshes,  when  compared  with  the  first  order  modified  high  order  schemes 
(the  right  panel  of  Figure  4.1). 
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5.4  Second  order  modification  to  the  high  order  discontinuous 
Galerkin  methods 

In  this  subsection,  the  second  order  entropic  projection  is  used  as  a  building  block  to  modify 
high  order  discontinuous  Galerkin  methods.  The  modification  follows  a  similar  line  as  the 
first  order  monotone  modification  presented  in  section  4.3.  The  procedure  for  a  (k  +  l)-th 
order  discontinuous  Galerkin  method  is  outlined  as  following. 

At  each  time  step  evolution,  we  would  like  to  update  Uj  =  {uf\  ■■■iuf\ulj,urj}  over  the 
cell  Ij.  At  the  initial  stage,  vtp  is  obtained  by  (2.10)  and  ul3  and  Uj  refer  to  the  point  values 
at  the  left  and  right  boundaries  of  the  cell  Ij.  The  scheme  consists  of  the  following  steps  in 
each  time  step  evolution. 

1.  Apply  the  WENO  limiters  in  the  oscillatory  troubled  cells  as  described  in  section  4.3. 

2.  Identify  nonconvex  troubled  cells  Ij  as  described  in  Criterion  II  in  section  4.3. 

3.  Identify  nonconvex  troubled  cell  boundaries  as  described  in  Criterion  I  in  section  4.2. 

4.  Modify  the  numerical  flux  ,f]+i  with  the  discontinuity  indicator. 

Specifically,  let  /•_,  i  =  t),  where  umf{  is  computed  from  (5.9)  at  a  troubled 

2  J~r  2  J~r  2  J~r  2 

cell  boundary.  At  good  cell  boundaries,  , . 

3'  2  3'  2 

5.  For  good  cells,  update  u(p  (/  =  0, ...,  k)  by  the  regular  discontinuous  Galerkin  method 
and  compute  ulj  and  urj  from  (l  —  0, ...,  k). 

6.  For  nonconvex  troubled  cells,  update  by  (4.5)  and  update  the  remaining  moments 
Uj'*  with  l  >  1  in  a  WENO  fashion  as  in  [12]. 

7.  For  nonconvex  troubled  cells  Ij,  update  ulj  and  urj  as  described  in  section  5.3. 

In  the  following,  we  provide  numerical  examples  of  the  modified  high  order  discontinuous 
Galerkin  methods  with  the  second  order  entropic  projection. 
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Example  5.7.  The  nonconvex  conservation  law  (4.4).  Table  5.3  gives  the  L 1  errors  and 
the  corresponding  orders  of  accuracy  for  the  modified  DG-WENO  scheme  with  the  entropic 
projection.  The  L 1  errors  of  the  modified  DG-WENO  method  are  comparable  to  those  of 
the  regular  DG-WENO  method  shown  in  Table  4.2. 


Table  5.3:  The  L 1  errors  and  the  corresponding  orders  of  accuracy  for  the  modified  DG- 
WENO  method  with  the  entropic  projection  at  the  time  t  =  0.2.  The  numerical  flux  is  the 
Godunov  flux. 


No.  of 

k  =  1 

k  =  2 

k  =  3 

points 

L 1  error 

order 

L 1  error 

order 

L 1  error 

order 

30 

5.313E-03 

1.146E-04 

6.951E-04 

60 

1.244E-03 

2.09 

4.924E-05 

4.54 

5.714E-06 

6.93 

90 

5.382E-04 

2.06 

1.264E-05 

3.35 

7.044E-07 

5.16 

120 

2.866E-04 

2.19 

5.120E-06 

3.14 

2.048E-07 

4.29 

150 

1.649E-04 

2.47 

2.608E-06 

3.02 

8.291E-08 

4.05 

Example  5.8.  The  nonconvex  conservation  law  (1.1)  with  f(u)  =  sin(u)  and  the  initial 
condition  (3.14).  The  solutions  at  the  time  t  —  4  are  plotted  in  Figure  5.3.  As  shown, 
the  numerical  solution  of  the  modified  high  order  discontinuous  Galerkin  method  with  the 
entropic  projection  successfully  converges  to  the  correct  entropy  solution. 

Example  5.9.  The  nonconvex  conservation  law  (1.1)  with  f(u)  defined  by  (3.5),  the 
initial  condition  (3.13)  and  a  periodic  boundary  condition.  As  shown  in  Figure  5.4,  the  high 
order  discontinuous  Galerkin  method  with  the  second  order  entropic  projection  approximates 
the  correct  entropy  solution  very  well  with  relatively  coarse  meshes.  When  comparing  with 
the  first  order  modified  high  order  discontinuous  Galerkin  method  (see  Figure  4.3),  the 
performance  of  the  scheme  has  been  greatly  improved. 

6  Concluding  Remarks 

We  have  investigated  the  performance  of  high  order  numerical  methods  for  general  one  di¬ 
mensional  scalar  nonconvex  conservation  laws,  emphasizing  convergence  to  the  discontinuous 
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Figure  5.3:  The  solid  line  is  the  exact  solution  of  (1.1)  with  (3.14)  at  the  time  t  —  4,  while 
the  solid  symbols  are  the  numerical  solution  of  the  modified  second  order  (left)  and  fourth 
order  (right)  discontinuous  Galerkin  method  with  the  entropic  projection.  The  uniform  mesh 
sizes  are  Ax  =  0.1  (unfilled  squares)  and  Ax  =  0.05  (solid  symbols),  respectively. 

entropy  solutions.  It  is  observed  that  high  order  finite  volume  WENO  and  discontinuous 
Galerkin  schemes  may  fail  to  converge  to  the  entropy  solutions  for  some  difficult  test  cases. 
A  first  order  modification  based  on  first  order  monotone  schemes  and  a  second  order  modifi¬ 
cation  based  on  an  entropic  projection  are  designed  for  high  order  finite  volume  WENO  and 
discontinuous  Galerkin  methods  to  maintain  high  order  accuracy  in  smooth  regions  and  to 
enforce  convergence  to  the  entropy  solution.  Numerical  examples  are  shown  to  demonstrate 
the  quality  of  the  proposed  schemes.  Future  work  will  include  more  efficient  implementa¬ 
tion  of  the  second  order  entropic  projection,  and  a  generalization  of  the  current  approach  to 
multi- dimensional  scalar  problems  and  to  systems. 
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